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The  discretisation  of  partial  differential  equations,  by  either  finite  element 
or  finite  difference  techniques,  often  leads  to  large  linear  systems  of  equa¬ 
tions  with  sparse  matrices.  Fast  iterative  solution  methods,  based  upon  the 
preconditioning  of  the  conjugate  gradients  method,  have  been  proposed  for  the 
symmetric  positive  definite  case and  also  for  more  general  situations.-Mi- 
vijy  In  this  report  we  present  new  sharp  upperbounds  for  the  conjugate  gradient 
residual.  These  upperbounds  help  us  to  understand  and  explain  the  convergence 
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behaviour  of  the  preconditioned  conjugate  gradients  method.  We  also  present  a 
type  of  preconditioning  that  has  almost  the  same  convergence  properties  as  those 
presented^  but  which  admit-  full  vectorization  on  supercomputers  like  the 

CRAY-1  and  the  CYBER  205.  ' 

For  the  nonsymmetric  case  we  propose  different  types  of  preconditioning  in 
connection  with  the  Chebyshev  iterative  method. 
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1.  INTRODUCTION 


The  aims  of  the  research,  as  carried  out  during  the  years  of  the  contracting 
period,  have  been  threefold: 

1.  To  explain  the  speed  of  convergence  of  the  preconditioned  conjugate  gradients 
method  (for  the  symmetric  positive  definite  case). 

2.  To  construct  fast  iterative  methods  for  sparse  linear  systems  with  nonsymmetric 
matrices. 

3.  To  develop  user-oriented  software  for  the  methods  mentioned  in  the  previous 
two  points. 

As  may  appear  from  section  2  and  3,  considerable  progress  has  been  made  with  respect 
to  all  the  above  listed  points.  The  results  which  have  been  obtained  have  all  been 
published  recently  or  they  have  been  submitted  for  publication.  Major  part  of  the 
research  has  also  been  published  in  the  form  of  a  Ph.D. -thesis  [3]. 

In  section  2  we  present  an  overview  of  the  results  obtained  for  the  synmetric 
positive  definite  case  and  in  section  3  we  introduce  the  results  for  the  nonsym¬ 
metric  case. 


2.  SYMMETRIC  POSITIVE  DEFINITE  MATRICES 

In  this  case  we  consider  the  conjugate  gradients  (eg)  method  for  the  iterative 
solution  of  K  ^Ax  *  K  *b,  where  Ax  ■  b  is  the  linear  system  to  be  solved  and  K  * 
is  a  preconditioning  matrix.  For  a  description  of  this  method  see  [  1  ) . 

In  (3,  chapter  2]  an  analysis  of  the  convergence  behaviour  of  the  eg  method  is 
presented,  which  exploits  its  relationship  to  the  Lanczos  method.  This  analysis 
leads  to  an  expression  for  the  residual  in  terms  of  the  Ritz  values  (i.e.  the 
eigenvalue  approximations  generated  by  the  Lanczos  method). 

It  has  been  observed,  e.g.,  by  Concus  et  al.  [4 1 ,  that  in  many  relevant 
situations  the  convergence  behaviour  appears  to  be  much  better  than  might  be 
expected  from  upperbounds  involving  only  the  condition  number  of  K  ^A.  This  so- 
called  superlinear  convergence  behaviour*  can,  even  for  rather  dense  spectra,  be 
fairly  well  explained  by  our  upperbounds  (see  (3,  chapter  2)),  which  take  advantage 
of  the  relative  separation  of  the  smallest  eigenvalues. 

For  a  specific  example  of  the  preconditioned  method  it  is  shown  that  the  sepa¬ 
ration  of  only  the  four  smallest  eigenvalues  is  sufficient  to  explain  the  actual 
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observed  convergence  behaviour.  Thi6  contrasts  with  the  widespread  belief  that  it 
is  the  clustering  of  the  eigenvalues  that  causes  the  fast  convergence. 

It  is  one  of  the  purposes  of  [3,  chapter  3]  to  show  that  for  a  class  of  model 
problems  (i.e.,  the  discretised  Poisson  equation)  the  preconditioning  of  A  should 
always  give  a  more  favourable  eigenvalue  distribution,  in  particular  with  respect 
to  the  smallest  eigenvalues. 

In  view  of  the  dependence  of  the  convergence  behaviour  of  eg  upon  the  separation 
of  the  smallest  eigenvalues  of  the  operator,  we  also  comnent  in  [3,  chapter  3]  on 
certain  modifications  of  the  ICCG  method,  which  have  been  proposed  in  litterature 
(e.g. ,  see  [5] ). 

In  [2]  a  family  of  incomplete  decompositions  is  presented,  to  be  used  as  pre¬ 
conditioning  matrices  for  the  eg  method.  The  question  which  member  of  this  family 
leads  to  optimal  results  (with  respect  to  efficiency)  remains  unanswered.  However, 
the  results  of  numerical  experiments  given  in  [2],  may  help  in  the  choice. 

For  the  simplest  types  of  incomplete  decomposition  Eisenstat  [6]  has  shown  that 
the  preconditioned  eg  algorithm  can  be  implemented  in  such  a  way  that  the  number 
of  arithmetical  operations  per  iteration  is  almost  the  same  as  for  the  unprecon¬ 
ditioned  eg  algorithm.  This  makes  the  use  of  these  decompositions  even  more  attrac¬ 
tive. 

The  use  of  incomplete  decompositions  on  vector  computers  gave  some  problems,  since 
they  do  not  automatically  lead  to  vectorizable  algorithms  1 7 , 8 ]  .  Only  with  the 
help  of  laborious  permutations  the  efficiency  of  this  algorithm,  as  compared  with 
other  algorithms,  can  be  restored  (8,9].  In  [10,11  ]  a  simple  variant  of  the  incom¬ 
plete  decomposition  algorithm  is  introduced  which  admits  full  vectorization  without 
notable  loss  in  the  number  of  iterations  and  instability. 

Software  for  the  preconditioned  eg  algorithm  is  available  in  the  program 
library  ACCULIB  [ 1 2 , 1 3  ]  . 


Remark.  For  the  actual  computation  of  the  extreme  eigenvalues  of  the  preconditioned 


3.  NONSYMMETRIC  MATRICES 


Here  we  consider  the  use  of  the  complex  Chebyshev  Iteration  method  for  the 
solution  of  K  'Ax  ■  K  *b,  where  K  '  is  a  suitable  preconditioning  matrix  for  the 
linear  system  Ax  ■  b.  This  method  converges  if  K  5 A  has  all  its  eigenvalues  in 
the  right  half  plane  [  1 5  ] . 

In  [16,18]  we  study  the  effects  of  preconditioning  on  the  Chebyshev  Iteration 
method,  when  applied  to  linear  sytems  that  arise  from  discretisation  by  finite 
differences  of  the  p.d.e. 
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over  a  rectangular  domain  in  R  ,  with  Dirichlet  boundary  conditions  and  with 
a(x,y)  >  0,  b(x,y)  >  0,  c(x,y)  >  0  . 

Different  preconditionings  K  *  are  considered,  most  of  them  are  based  upon 
incomplete  decomposition  of  A.  For  some  typical  examples  it  is  shown  that  the  latter 
type  of  preconditioning  leads  to  a  highly  competitive  iterative  method  as  compared 
to  a  number  of  other  methods.  For  many  relevant  situations  we  are  able  to  prove 
that  K  *A  has  all  its  eigenvalues  in  the  right  half  plane. 

Special  problems  may  be  encountered  when  the  coefficients  d(x,y)  or  e(x,y) 
in  the  p.d.e.  assume  large  values.  It  then  often  appears  that  the  factors  L  and  D 
of  the  incomplete  decomposition  are  very  ill-conditioned,  whereas  A  is  fairly  well- 
conditioned.  In  [16,17,3:  chapter  4]  we  show  how  the  ill-conditioning  of  the  factors 
can  be  remedied  by  certain  modifications  of  the  incomplete  decomposition  algorithm. 
These  modifications  involve  either  a  single  parameter  o  or  simply  replace  the 
diagonal  elements  of  L  and  U  by  values  that  guarantee  well-conditioning. 

We  can  prove  that  for  a  class  of  problems  there  are  values  of  o0  (depending  on 

-1  .  u 
the  operator),  such  that  K  A  has  all  its  eigenvalues  in  the  right  half  plane  for 

0  >  Cq.  An  unsolved  problem  remains  the  experimental  observation  that  the  smallest 

value  of  o  which  gives  well-conditioned  factors  is  also  optimal  with  respect  to 

the  6peed  of  convergence. 

Analogously  to  the  preconditioned  eg  method  (see  section  2)  the  complex  Chebyshev 
Iteration  method  can  be  formulated  in  such  a  way  that  for  the  most  relevant  precon¬ 
ditionings,  the  preconditioned  iterations  require  only  little  more  computer  time 
than  the  unpreconditioned  iterations  (see  [)6]). 

Variants  of  the  incomplete  decomposition  preconditionings  that  admit  full  vecto- 
rization  on  vector /parallel  processors,  can  be  easily  constructed  as  is  shown  in 
[10,11]. 


Software  for  the  preconditioned  complex  Chebyshev  Iteration  method  is  avai¬ 
lable  in  the  program  library  ACCULIB  [13]. 
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